#这是工程实例计算(分段二次牛顿插值多项式拟合单变量工程用表数据并自动化)
import numpy as np
import pandas as pd
from matplotlib import pyplot as plt
import sys

class n2Interpolate_multinomial():
    def __init__(self, path, x, path_s, **kwargs):
        """
        这是牛顿法分段二次插值多项式。
        :param path: 文件路径，csv格式(x,y),数据对数量最好是大于等于3的奇数。
        :param x: 计算的x值
        :param kwargs: 坐标轴名称 xlabel=“”，ylabel=“”, tittle=“”，拟合步长step=,
        数据理论最大值最小值max=,min=,默认所给数据的最大最小值,
        :param path_s: 保存路径path=""
        """
        self.pa = False #打印结果标志
        self.x = x
        self.path_s = path_s
        try:
            self.path = path
            self.step = kwargs["step"]
        except:
            print("参数输入不全！")
            sys.exit(1)
        try:
            self.t = kwargs["tittle"]
        except:
            self.t = ""
        self.data = pd.read_csv(path)
        if self.data.shape[0] % 2 == 0:
            self.data = self.data.drop(axis=0, index=self.data.shape[0]-1)
        self.data = np.array(self.data)
        self.f = []
        self.ch_data = []
        for i in range(0, self.data.shape[0], 2):
            if i != self.data.shape[0] - 1:
                self.ch_data.append(self.data[i:i+3])
        try:
            self.max = kwargs["max"]
            self.min = kwargs["min"]
        except:
            self.max = self.data[:, 1].max()
            self.min = self.data[:, 1].min()
        self.caculate(self.x)
        try:
            self.paint(x=kwargs["xlabel"], y=kwargs["ylabel"], step=kwargs["step"], tittle=self.t)
            self.save()
        except:
            print("参数输入有误！")
            sys.exit(1)

    def forward(self, data):
        self.difference_quotient = np.zeros((data.shape[0], data.shape[0] + 1))
        for j in range(data.shape[0]):
            self.difference_quotient[j][0] = data[j][0]
            self.difference_quotient[j][1] = data[j][1]
        for i in range(2, data.shape[0] + 1):#列
            for j in range(i-1, data.shape[0]):#行
                self.difference_quotient[j][i] = (self.difference_quotient[j][i-1] - self.difference_quotient[j-1][i-1])/(self.difference_quotient[j][0] - self.difference_quotient[j-(i-1)][0])
        self.fdq = np.delete(self.difference_quotient, 0, axis=1)
        return self.fdq

    def caculate(self, x):
        hanshu = []
        result = []
        for i in range(len(self.ch_data)):
            hanshu.append(self.forward(self.ch_data[i]))
        for i in range(len(hanshu)):
            self.y = hanshu[i][0][0] + hanshu[i][1][1]*(x - self.ch_data[i][0][0]) + hanshu[i][2][2]*(x - self.ch_data[i][0][0])*(x-self.ch_data[i][1][0])
            result.append(self.y)
        self.result = result
        flag = False
        for i in range(len(self.ch_data)):
            if x<self.ch_data[i][2][0] and x>=self.ch_data[i][0][0]:
                self.result = self.result[i]
                flag = True
            elif x<self.ch_data[0][0][0] and x>self.ch_data[-1][2][0]:
                print("计算值超限！")
                sys.exit(1)
            elif x == self.ch_data[-1][2][0]:
                self.result = self.ch_data[-1][2][1]
                flag = True
        if self.pa == False:
            if flag == True:
                print("x={}时，结果是{}。".format(x, self.result))
                return self.result
            else:
                print("计算值超限！")
                sys.exit(1)
        else:
            if flag == True:
                return self.result
            else:
                print("计算值超限！")
                sys.exit(1)

    def paint(self, x, y, step, tittle=""):
        self.pa = True
        self.fig = plt.figure(1, figsize=(10, 5))
        plt.rcParams["font.sans-serif"] = ["SimHei"]
        plt.rcParams["axes.unicode_minus"] = False
        plt.scatter(self.data[:, 0], self.data[:, 1], marker="o", c="red")
        plt.grid()
        plt.xlabel(x, fontsize=15)
        plt.ylabel(y, fontsize=15)
        plt.title(tittle, fontsize=15)
        xx = []
        yy = []
        a = self.data[0][0]
        b = self.data[-1][0]
        self.xx = np.arange(a-0.1, b, step)
        for i in range(self.xx.shape[0]):
            if self.xx[i] <= self.data[-1][0] and self.xx[i] >= self.data[0][0]:
                self.x = self.xx[i]
                xx.append(self.xx[i])
                if self.caculate(self.x) > self.max:
                    yy.append(self.max)
                elif self.caculate(self.x) < self.min:
                    yy.append(self.min)
                else:
                    yy.append(self.caculate(self.x))
            else:
                pass
        plt.plot(xx, yy, c="blue")
        plt.show()

    def save(self):
        self.fig.savefig(self.path_s + ".png")

if __name__ == '__main__':
    #饱和水饱和蒸汽
    n2Interpolate_multinomial("./data/饱和蒸汽表.csv", 1.44, xlabel="相对气压（巴）", ylabel="饱和温度（℃）", step=0.01,
                              path_s="./data/饱和蒸汽表牛顿分段插值", tittle="压强-饱和温度")

    # #单位：m，mH20,海拔插值
    # n2Interpolate_multinomial("./data/海拔压力表.csv", 0, xlabel="海拔（m）", ylabel="压强（mH20）", step=10, path_s="./data/海拔压力牛顿分段插值", tittle="压强-海拔图")
    #
    # #标准正态分布表 P(X<=x)=y
    # Normal_distribution = n2Interpolate_multinomial("./data/zmj.csv", 1.88, xlabel="x", ylabel="y", step=0.01, path_s="./data/zmj_save", tittle="标准正态分布")
